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1 Introduction 

The recent biological theory of adaptive dynamics [6,9] proposes a descrip- 
tion of the long term evolution of an asexual population by putting emphasis 
on the ecological interactions between individuals, in contrast with classi- 
cal population genetics models which focus on the genetic structure of the 
population. The basic models are individual-based models in which the 
population dynamics is precisely described and includes birth, death, com- 
petition and mutation. The basic idea of the theory of adaptive dynamics 
is to try to get insights about the interplay between ecology and evolution 
by studying the invasion of a single mutant type appearing in a simplified 
(monotype stable) resident population. The evolution of the population can 
then be described as a sequence of mutant invasions in the population. If 
the resident type goes extinct when the mutant type invades (we say that 
the mutant type fixates), the evolution is described by the so-called 'trait 
substitution sequence' (TSS) [10]. This appraoch has revealed powerful to 
predict the qualitative behaviour of complicated evolutionary dynamics. In 
particular, it allows to determine the (local) direction of evolution in the 
space of phenotypic traits (or simply traits) from the individual ecological 
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parameters, and to predict and explain the phenomenon of evolutionary 
branching [10], where a population, initially composed of individuals with 
traits concentrated around a single trait value, divides into two (or more) 
subpopulations concentrated around distinct trait values that stably coexist 
because of their ecological interactions. The description of this phenomenon 
is an important achievement of this theory, as well as the 'canonical equa- 
tion of adaptive dynamics' [3], which describes the evolution of the dominant 
trait of the population as a deterministic "hill-climbing" process on a fitness 
landscape which depends on the current state of the population (see ([I]) 
below) . 

More formally, as soon as eternal coexistence of two or more types is not 
permitted by the model, the evolution can be described by the sequences 
(Tn)n and {Vn)n, where is the n-th time where the population becomes 
monomorphic (i.e. composed of only one type) and Vn is the surviving type 
at time T„. The sequence (V^) is the above-mentionned TSS. It is possible 
to prove the convergence of an individual-based model to the TSS under 
two biologically motivated assumptions [10,1]. First, the assumption of rare 
mutations guarantees that, in the timescale of mutations, the widths of 
time intervals during which the population is polymorphic vanish, so that 
there is only one type surviving at any time t. To prevent the population 
from rapidly becoming extinct in the new timescale, one also has to rescale 
population sizes, thereby making the assumption of large populations. 

Subsequently, the TSS is a Markov jump process on the trait space whose 
semigroup is shown [1] to depend on the invasion fitnesses (as defined in [9]) 
f{x,y), x,y G X, where f{x,y) is defined as the expected growth rate of 
a single individual of type y — the mutant — entering a monomorphic 
population of type x 'at equilibrium' — the residents. Note that this fitness 
is not given a priori, but derived from the microscopic model of individual 
interactions. Because of the assumption of large population, the sign of this 
fitness determines the possibility of invasion of a mutant type: if f{x, y) < 0, 
the mutant type y cannot invade a resident population of type x. Thus, 
evolution proceeds by successive invasions of (only) advantageous mutant 
types replacing the resident one. 

The 'canonical equation of adaptive dynamics' [3], which describes the 
evolution of a one-dimensional trait x as the solution of the following ODE, 
is obtained from the TSS in the limit of small mutations: 





where 




stands for the (rescaled) variance of the mutation step law. 
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n{x) for the equihbrium size of a pure x-type population, and f(x,y) for the 
invasion fitness mentioned above. Note how only advantageous types get 
fixed (the trait follows the fitness gradient) and how the fitness landscape 
y !—>■ f(x, y) depends on the current state x of the population. 

However, it is well-known that slightly deleterious types can be fixed by 
chance in finite populations. This phenomenon is known under the name of 
genetic drift. Depending on the strength of genetic drift, selection is said to 
be strong (genetic drift has negligible effects) or weak. In the large popula- 
tion asymptotic from which the TSS of adaptive dynamics is derived, genetic 
drift has negligible impact compared to the action of selection. Therefore, 
the fixation of slightly deleterious types cannot be observed. Our goal here 
is to include genetic drift in the adaptive dynamics models by considering 
finite populations under weak selection. We continue using the bottom-up 
approach of adaptive dynamics; that is, model (macroscopic) evolution from 
(microscopic) populations. In particular, we allow the population sizes to 
fluctuate randomly through time and we aim to reconstruct a fitness function 
from the microscopic parameters. 

After the description of the model (Section 2), we derive a new TSS in 
the limit of rare mutations (Section 3), from which a limit of small mu- 
tations gives what we call the 'canonical diffusion of adaptive dynamics' 
(Section 4). The coefficients of this diffusion involve the first-order deriva- 
tives of the fixation probabilities, which are computed in Section 5 as a linear 
combination of four fundamental components associated to fertility, defence, 
aggressiveness and isolation. New numerical results on the robustness of the 
population with respect to these fundamental components are also given, as 
well as sonic consequences on the canonical diffusion of adaptive dynamics 
in large populations. 

2 The microscopic model 

We will restrict here to logistic interaction. More general models are con- 
sidered in [2]. 

A monotype (binary) logistic branching process (LBP, see [7]) with dy- 
namical parameters (6, c, d) is a Markov chain in continuous time {Xf, t > 0) 
with nonnegative integer values and transition rates 




bi if 
ci{i — 1) if 
-i{b + c{i-l)) if 







otherwise. 



j = i + l 
j = i-l 
j = i 
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The non-linear term ci{i — l) describes competition mortality due to random 

encounters between individuals. Other terms correspond to independent 
birth events with constant individual rates. This Markov chain is positive- 
recurrent and converges in distribution to a r.v. ^, where ^ is a Poisson 
variable of parameter 9 := b/c conditioned on being nonzero 

e"^ 9'^ 

m = ^) = Y^-( ^>i. (2) 

We consider in the sequel, a multitype asexual birth and death pro- 
cess with mutation, generalizing this LBP. At any time t, the population is 
composed of a finite number N{t) of individuals characterized by their phe- 
notypic traits xi{t), . . . , Xjv(t) (t) belonging to a given trait space X, assumed 
to be a closed subset of MJ^. The population state at time t is represented 
by the counting measure on X 

N{t) 
i=l 

The population dynamics is governed by the following parameters. 

• b{x) is the rate of birth from an individual of type x. The function b 

is assumed to be C^. 

• c(x, y) is the rate of death of an individual of type x due to the com- 
petition with another individual of type y. Therefore, the total death 
rate of an individual of type a; in a population u may be written as 

J c{x,y){v{dy) — 5x{dy)). In this expression, the Dirac mass at x sub- 
stracted to the measure u means that the individual docs not compete 
with himself. The function c is assumed to be C| and bounded away 
from on X'^. 

• ■~ilJi{x) is the probability that a birth from an individual with trait x 
produces a mutant individual, where /Lt(x) G [0, 1] and where 7 G (0, 1] 
is a parameter scaling the frequence of mutations. When there is no 
mutation, the new individual inherits the trait of its progenitor. In 
Section 3, we will be interested in the limit of rare mutations (7^0). 

• M(x, dh) is the law of the trait difference h = y — x between a mutant 
individual with trait y born from an individual with trait x. We assume 
that M{x, dh) has expectation (no mutation bias), i.e. / hM{x, dh) = 
0, and has a density on which is uniformly bounded in x € by 
some function M{h) with finite third-order moment. 
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We will denote the dependence of vt on the parameter 7 with the notation 
. Observe that such a population cannot go extinct because the death 
rate is when there is only one individual in the population. Since we want 
to apply a limit of rare mutations while keeping the population size finite, 
this is necessary to prevent the population to become extinct before any 
mutation occur. 

Let ^(x) be a random variable whose law is the stationary distribution 
of a pure x-type population with no mutation = 0). This law is given 
by ([2]) where is replaced by 0{x) := b{x)/c{x,x). 

The last notation needed concerns a population with initially only two 
types X and y and with no mutation. Then vt = ^t^x + Yt^y, where {Xt, Yt : 
t > 0) is a bivariate Markov chain. For this Markov chain, P(r < 00) = 1, 
where T is the first time where either Xt or Yt reach 0. We call fixation 
(of the mutant y) the event {Xt = 0}. The probability of fixation will be 
denoted by Un,mix,y) 

Un,mix,y) ■■= F{Xt = I Xo = n,yo = in). 

3 The trait substitution sequence in finite popu- 
lations 

In this section, we apply the limit of rare mutations (7 0) to the process 
I'"', in order to describe the evolution of the population as a TSS in finite 
population. This limit requires to rescale time properly, as t/^, to describe 
the evolution on the mutation timescale. 

Theorem 3.1 Fixx G X. Assume that i'q = Nq6x where sup ^^(^q^i-^M{{NqY) < 
00 for somep > 1. Then, for any < ti < . . . < the n-tuple (yj^j^^ . . . , 
converges in law for the weak topology to {Nt-^5st^ ; • • • j ^t„Sst„) where 

(1) {St;t >0) is a Markov jump process on X with initial value Sq = x and 
whose jumping rates q{x, dh) from x to x + h are given by 

q{x, dh) = I3{x)xix, x + h)M{x, dh), 

where /?(x) = ^(x)6(x)E(^(x)) = ^(x)6(x)6'(x)/(l - e"^^""^) and 

, . ^nP(e(x)=n) ^ , n 

(3) 
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(2) Conditional on {Sf-^, . . . , St„) = {xi, . . . ,Xn), the are independent 
and respectively distributed as (,{xi). 

Therefore, in the hmit of rare mutations, on the mutation timescale, the 
population is always monomorphic and the dominant trait of the population 
evolves as a jump process over the trait space, where a jump corresponds to 
the appearance and fixation of a mutant type. Moreover, at any time, the 
population size stationary (i.e. has the stationary distribution corresponding 
to the dominant trait of the population). The fixation rate of a mutant is 
governed by the function xi^-, y)i which is therefore the random analogue of 
the traditional invasion fitness [9], defined as the probability of invasion of 
a mutant type y is a resident population of type x at equilibrium. Observe 
that, as usual in adaptive dynamics, the fitness landscape depends on the 
current state of the population. Moreover, in contrast with the classical 
TSS [10], from a given monomorphic resident population, any mutant trait 
has a positive probability to invade (by genetic drift). Therefore, evolution is 
possible in any direction of the trait space. However, a directional selection 
still exists, as will appear in the next section. 

We refer to [2] for the proof of this result. However, this convergence is 
natural in view of the following interpretation of each parameter. f3{x) can 
be seen as the mean mutant production rate of a stationary x-type population 
(i.e. with size ^(.x)), and x(a;,y) is the probability of fixation of a single y- 
type mutant entering a pure x-type population with size-biased stationary 
size. The size bias comes from the fact that the mutant appears at a birth 
time in the stationary population (since the birth rate is proportional to 
the population size, the population size after a birth event in the stationary 
population is given by the size-biaised stationary population size). 

4 The canonical diffusion of adaptive dynamics 

Let us assume for simplicity that ^ = M^. Let a{x) be the square root 
matrix of the covariance matrix of M{x, •). We also need to assume that the 
matrix ^(x) is a Lipschitz function of x. 

In order to obtain the equivalent of the canonical equation of adaptive 
dynamics in a finite population, we want to apply a limit of small mutation 
steps (weak selection) to the TSS S. To this aim, we introduce a param- 
eter e > and replace the mutation kernels M(x, •) by their image by the 
application h eh. Time also has to be rescaled in order to obtain a 
non-degenerate limit. The correct time scaling is which leads to the 
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following generator for the rescaled TSS {S^; t > 0) 

A^(f{x) = -^ I {(f{x + eh) — (f{x))P{x)x{x,x + eh)M{x,dh). 

Using the assumption that the mutation kernels M(x, •) have expectation, 
it is elementary to compute the limit of this expression as e — > (for suffi- 
ciently regular ip). This limit, which takes the form of a diffusion generator, 
explains the following result (its full proof can be found in [2]). 

Theorem 4.1 // the family (5'Q)e>o has bounded first-order moments and 
converges in law as e ^ to a random variable Zq, then the process 5^ with 
initial state Sq converges in law for the Skorohod topology on D(IR_(.,R^) to 
the diffusion process {Zt;t > 0) with initial state Zq unique solution to the 
stochastic differential equation 

dZt = (3iZt)a\Zt) ■ V2X{Zu Zt)dt + ^ (3{Zt)x{Zu Zt)a{Zt) ■ dBt (4) 

where V2X denotes the gradient w.r.t. the second variable y ofx{x,y) and 
B is a standard k- dimensional Brownian motion. 

This result gives the equivalent of the canonical equation of adaptive 
dynamics ([T]) when the population is finite. It is no longer a deterministic 
ODE, but a diffusion process, in which the genetic drift remains present (in 
the form of a stochastic diffusion term), as a consequence of the population 
finiteness and of the limit of weak selection. The deterministic drift part 
of ^ is very similar to the standard canonical equation of adaptive dynam- 
ics ([1]), and involves in particular the gradient of the fitness function x- The 
process ([4]) provides a diffusion model describing the evolution of the dom- 
inant trait value in a population [8,5], grounded on a precise microscopic 
density-dependent modelling of the population dynamics. It also gives the 
precise balance between directional selection and genetic drift as a function 
of the individual's dynamical parameters. 

5 Fixation probability near neutrality 

The SDE (j4j) involves the fixation probability x{x-,x) and the fitness gradient 
with respect to the second variable V2X{x-,x). In this section, we explain 
how these quantities can be explicitly computed. 

We need to compute the derivatives of the fixation probabilities Un,m{x^ y) 
when y is close to x. Recall that the law of the two- types LBP without mu- 
tation (X, Y) used to define Un,m in the end of Section 2 is characterized by 
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the birth vector B and the competition matrix C 



B={ C-^ c{x,x) c{x,y) 



b{y) J ' V c{y,y) 

We will say that the mutant is neutral if all individuals are exchangeable, 
i.e. when b{y) = b{x) and c(x, y) = c{y, x) = c(y, y) = c{x, x) (this holds in 
particular when y = x). As will appear below, using the notation (6, c) := 
{b{x) , c{x , x)) , it is natural to focus on deviations from the neutral case 
expressed as 

B = bl+('), C = cl-^' ^0 ^0 ' 



Ay' \ 6 S J \ a J \ £ 

In words, deviations from the neutral case are a linear combination of four 
fundamental selection coefficients A, S, a, e, that arc chosen to be positive 
when they confer an advantage to the mutant. It is convenient to assess 
these deviations to the neutral case in terms of 

1. fertility (A): positive A means increased mutant birth rate 

2. defence capacity (6): positive S means reduced competition sensitiv- 
ity of mutant individuals w.r.t. the total population size 

3. aggressiveness (a): positive a means raised competition pressure 
exerted from any mutant individual onto the rest of the population 

4. isolation (e): positive e means lighter cross-competition between the 
two different types, that would lead, if harsher, to a greater probability 
of exclusion of the less abundant one 

Under neutrality, an elementary martingale argument shows that the fixa- 
tion probability equals the initial mutant frequency p := m/{m + n). This 
implies in particular that 

e-^W-l + g(x) 
X{x,x) = . (5) 

The following theorem unveils the dependence of u upon A, 5, a, e, when 
they slightly deviate from 0, and explains why these four selection coefHcients 
provide a natural basis to decompose the gradient of the fixation probability. 
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Theorem 5.1 As a function of the multidimensional selection coefficient 
s = {X,6,a,£), the probability u is differentiable, and in a neighbourhood of 
s = (selective neutrality), 



u = p + v' .s + o{s) 



(6) 



where the selection gradient v = {v^,v ,v'^,v'^) can be expressed as 



n,m 



n,m 



P{'^-P) 9n+m 
p{l-p){l- 2p) 5^+^ 



i ^ 8, 



where the g's depend solely on the resident's characteristics b,c, and on the 
total initial population size n+m. They are called the invasibility coefficients. 

The invasibihty coefficients of a pure resident population are interesting 
to study, as they provide insights about how the fixation probabihty devi- 
ates from p as the selection coefficients of the mutant deviate from 0. They 
also provide information about the robustness of the resident population, 
i.e. its resistance to mutant invasions. In particular, this allows to compare 
the sensitivity of the invasion probability in a given monomorphic resident 
population with respect to the four fundamental selection coefficients. In 
the simplest case where mutations in the parameter space are isotropic, 
the biggest invasibility coefficent gives the direction of the parameter space 
where a mutant is more likely to invade. More generally, when there are 
correlations between mutations in the parameter space (either because of 
the phenotypic structure in the functions b[-) and c(-, •), or because the co- 
variance matrix of the mutation steps is non-diagonal) , the likeliest direction 
of evolution in the trait space is given by the deterministic coefficient of the 
canonical diffusion ([3]), in which the fitness gradient is given by 

V2x{x, x) = a\{x)Vb{x) - as{x)'\/ ic{x , x) + aa(x)V2c(x, x), (7) 

where, for i = A, (5, a. 



It is possible to obtain explicit expressions for the invasibility coefficients 
g'' as series. We refer to [2] for the exact expressions. In particular, these 
expressions yield that at(x) = a^{6{x))/c{x, x) for some function a^. More- 
over, they allow one to compute numerically the invasibility coefficients, and 
therefore the quantities for l = \,5,a as functions of the parameter 9{x). 
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Figure 1: The functions ax, as and functions of 0. 




Figure 2: The functions 9 i— > 6a\{6) and 6 i— > 6aa{0). 
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These numerical results can be used to make simulations of the canonical 
diffusion of adaptive dynamics in various ecological examples. In particular, 
in contrast with the classical canonical equation of adaptive dynamics, the 
presence of a genetic drift can induce the evolutionary dynamics to drift 
away from evolutionary stable strategies, where the fitness gradient is zero. 
When the fitness gradient admits several zeros, it can visit several basins of 
attraction on various timescales. 

This numerical study is a work in progress, that is quite delicate because 
the series involved in the computation of g'' are slowly converging, with first 
terms that grow exponentially fast with 9. We shall give here our first 
results. Fig. [T] shows the functions for l = \,5,a. Several comments can 
be made from this figure. First, for any 9 > 0, as{9) > da{9) > dx{9). This 
means that, for equal mutation steps in the parameter space, a mutation is 
always more advantageous in the direction 5 than in the direction a, which 
is itself more advantageous than in the direction A. In other words, in a 
given population, a better defence capacity is more beneficial than a better 
aggressiveness, which is more beneficial than a better fertility. 

Moreover, as 9 goes to infinity, these functions have different asymptotic 
behaviors, ds seems to converges to 1/2, whereas dx{6) and da{9) are both 
equivalent to 1/29 (see Fig. [2]). However, this does not mean that mutations 
are much more likely to fixate in the 6 direction, because, when 9 is large, 
b = 9c is larger than c, and, in ([7]), and Oq, are multiplied by Vic and V2C 
respectively, whereas ax is multiplied by V6. More formally, to compute 
the limit of the canonical diffusion when 9 goes to infinity, one can divide 
the competition kernel c(-, •) by a constant K in the microscopic model and 
then let K go to infinity in the canonical diffusion. Denoting by xxi', ■) the 
fitness function obtained this way and using the asymptotic behaviors given 
above, one gets that 



Moreover, by ([5]), XKix,x) converges to when K ^ 00. Now, as proved 
in [1], with our notation, the fitness function of the canonical equation of 
adaptive dynamics ([T]) is given by f{x,y) = b(y) — c{y, x)9{x). Therefore, 
as K ^ +00, the canonical diffusion converges to a deterministic ODE 
which is precisely the canonical equation of adaptive dynamics. This gives 
a new justification of this equation, and this also allows one to study the 
fluctuations around the canonical diffusion when K is large. In particular, 
this diffusion equation with small diffusion term enters the framework of 
Freidlin-Wentzell's theory [4], which can be applied to predict the long time 



lim W2Xk{x,x) 



1 



{Vb{x) -9{x)Vic{x,x)). 
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behaviour of the diffusion, and its chain of visit of basins of attractions when 
K is large. This kind of information is biologically very relevant, since it 
allows one to predict in which order all the evolutionary stable strategies 
will be visited by the population and on which timescale. 
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